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AN  ANALYSIS  OF  THE  FLOW  FIELD  AROUND  A  2-D 
BODY  OF  ARBITRARY  SHAPE 

« 

I .  I ntroduction 

Purpose 

The  purpose  of  this  study  was  to  determine  the  velocity 
field  and,  therefore,  the  pressure  distribution  about  a  body 
of  arbitrary  shape  in  a  two-dimensional,  steady,  incompress¬ 
ible,  constant  viscosity  flow,  utilizing  a  combination  of 
Theodorsen's  method  and  the  Karraan-Pohlhausen  method.  A 
computer  program  to  implement  this  combined  method  has  been 
written,  and  results  from  this  program  are  described  herein. 

Some  practical  and  theoretical  considerations  leading 
to  this  combined  method  are  presented  in  the  remainder  of 
this  section. 

Background 

The  theoretical  analysis  or  design  of  an  airfoil  shape 
is  based  primarily  upon  a  knowledge  of  the  velocity  field 
in  the  neighborhood  of  the  body.  Since  the  time  of  Prandtl, 
aerodynamic! sts  have  realized  that  the  flow  over  a  body 
could  be  divided,  for  purposes  of  analysis,  into  two  regions 
Observations  indicate  that  viscous  forces  predominate  over 
inertia  forces  in  a  very  thin  layer  near  the  surface  of  the 
body,  and  this  region  has  come  to  be  called  the  boundary 
layer.  Outside  this  region,  the  flow  may,  with  little 
error,  be  considered  inviscid.  Since  inviscid  flow  with  a 
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uniform  parallel  onset  velocity  may  be  described  in  terms 
cf  a  velocity  potential,  the  inviscid  region  is  also  called 
the  region  of  potential  flow.  If  the  velocity  field  can 
be  determined  in  the  inviscid  flow,  then  the  velocity  at 
the  edge  of  the  boundary  layer  will  be  known  and  may  be 
used  as  a  boundary  condition  for  solution  of  the  boundary 
layer  equations.  Each  of  these  regions  will  now  be 
considered . 

Region  of  Inviscid  Flow 

For  thin  streamlined  airfoils  at  low  angles  of  attack, 
small  perturbation  theory  permits  a  number  of  simplifying 
assumptions  in  the  equations  of  motion  for  inviscid  flow. 
However,  no  such  simplified  potential  theory  was  available 
for  airfoils  of  arbitrary  shape  until  1933,  when  Theodorsen's 
method  (Ref  2)  was  published. 

Theodorsen  made  use  of  the  fact  that  any  closed  curve, 
-"'“h  as  an  airfoil,  defined  in  the  complex  plane,  may  be 
transformed  mathematical ly  into  a  circle.  This  transforma¬ 
tion  may  be  further  required  to  be  conformal;  that  is,  to 
preserve  the  local  angular  relationship  between  lines  passing 
through  each  point.  Conformality,  therefore,  ensures  that 
the  angle  of  attack  of  the  airflow  is  unaffected  by  the 
transformation.  The  problem  then  becomes  that  of  analyzing 
the  two-dimensional  flow  over  a  circle.  Since  this  is  one 
of  the  very  few  cases  for  which  exact  solutions  to  the 
equations  of  motion  are  known,  the  velocity  at  each  point 
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on  the  circle  nay  be  obtained  immediately.  Then,  through 
an  inverse  transformation,  the  potential  flow  velocity  at 
the ; corresponding  point  on  the  surface  of  the  airfoil  may 
also  be  determined.  This  velocity  may  be  used  as  an 
approximation,  (since  boundary  layer  thickness  has  not  yet 
been  accounted  for)  to  the  velocity  at  the  edge  of  the 
boundary  layer. 

Straightforward  though  this  method  appears,  its  use 
before  the  advent  of  high-speed  computers  was  linited  in 
practice  by  the  enormity  of  the  task  of  solving  the  necessary 
integral  equations  iteratively  to  the  desired  accuracy. 

The  transformation  is  very  tiae-consuming  and  tedious  when 
performed  by  hand,  providing  numerous  opportunities  to  make 
minor  mathematical  errors  which  may  invalidate  subsequent 
work.  Since  these  are  precisely  the  disadvantages  in 
computation  that  a  digital  computer  is  designed  to  overcome, 
Theodorsen's  method  lends  itself  well  to  formulation  for 
machine  use. 

The  rate  at  which  the  iterative  solution  converges 
is  largely  dependent  upon  how  nearly  the  original  airfoil 
approximates  a  circle.  Therefore,  Theodorsen  recommends 
th~  use  of  an  intermediate  analytical  transformation  (the 
Joukowski  transformation)  which  converts  the  airfoil  shape 
into  a  pscudo-circlc  conformally,  and,  thereby,  reduces 
the  computation  tine  required. 
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Boundary  Layer  Reg ?  <v« 

The  boundary  layer  equations  nay  be  derived  from  the 
Navier-Stokes  equations  for  viscous  fluid  notion  under  the 
assumptions  of  steady,  incompressible,  constant  viscosity 
flovr  at  large  Reynold's  numbers.  The  Karaan-Pohlhausen 
method  (Ref  3)  of  solving  the  boundary  layer  equations  is 
based  on  the  further  assumption  that  the  velocity  profile 
may  be  adequately  represented  as  a  fourth-degree  polynomial 
satisfying  appropriate  boundary  conditions,  including  the 
edge  velocity  distribution  previously  found  from  Theodorsen's 
method . 

Combined  F 1  ov;  Fields 

The  boundary  layer  flow  determined  from  the  Karman- 
Pohlhausen  method  is  only  a  first  approximation,  because  a 
body  with  a  boundary  layer  appears  to  the  potential  flow  to 
be  thicker  by  the  amount  of  the  displacement  thickness. 
Therefore,  the  displacement  thickness  of  the  boundary  layer 
must  be  computed  for  the  approximate  edge  velocity;  then  this 
thickness  must  be  added  to  the  actual  thickness  of  the  body; 
and  the  inviscid  velocity  field  for  the  resulting  shape 
determined  by  Theodorsen's  method.  Thus,  the  velocity  field 
solutions  are  continually  refined  until  the  change  in  dis¬ 
placement  thickness  from  one  iteration  to  the  next  is 
negligible.  At  this  point,  the  problem  of  finding  the 
velocity  field  nay  be  considered  solved.  Viith  the  velocity 
field  known,  other  aerodynamic  parameters  such  as  lift, 
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drag,  circulation,  and  t  .'par  at  ion  point  eay  also  be  computed, 
using  standard  sethods. 

1 
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1 I .  Thcodorsen's  Method 

Method  of  Presentation 

Theodorsen's  sethod  of  irrotational  flow  analysis 
depends  upon  several  inportant  results  of  potential  theory 
and  the  theory  of  conformal  transformations  of  coupler, 
functions.  These  results  Kill  be  briefly  sunaarized  before 
a  detailed  presentation  cf  the  Joukouski  and  Theodorsen 
transforaations  is  given.  Finally,  the  process  for  obtaining 
surface  velocities  by  an  inverse  transformation  will  be 
presented  - 

Review  of  Potential  Theory 

The  flow  outside  the  boundary  layer  is  ass used  to  be 
steady,  irrotational,  ar.d  inconpressible ,  permitting  its 
description  in  terns  of  a  velocity  potential  /  or  a  stress 
function  defined  so  that: 

..  _  30  _  3  r  .  ,  -  . 

u  -  7—  *  velocity  m  x-direction 

cx  3y 

(1) 

S*  3  V 

v  =  jy  =  -  53^  velocity  in  y-direction 

The  lines  of  constant  V  (streanlines)  and  lines  of  constant 

4  (equipotent ia 1  lines)  are  orthogonal  families  of  curves. 

Any  functions  of  C>  and  Y  which  satisfy  Laplace's 
a  •> 

equation,  =  V“Y  =  0,  represent  possible  types  of  fluid 

notion.  The  boundary  conditions  of  i mnernonbi  1  i  ty  of  lb«* 
airfoil  surface  and  uniform  flow  at  infinity  define  a  unique 
solution. 
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Con  f  or»  a  1  Trans foraatioa 

If  v(()  s  ♦(x,y)  ♦  i?(x,y)  is  aa  analytic  foactioa  of 
tke  coupler  variable  5  *  x  ♦  iy,  and  the  coapoaeats  i  and  T 
satisfy  Laplace*'  equation,  then  certain  useful  relations 
exist  between  the  C-plane  and  the  w-plane.  Specifically, 
any  sinple  curve  f(£)*  such  as  an  airfoil  surface,  naps  into 
a  curve  f(w)  ia  the  w-plane  in  such  a  nanner  that  the  angles 
between  pairs  of  lines  passing  through  any  point  in  the 
C-plane  renain  unchanged  at  the  inage  point  in  the  w-plane, 
although  local  rotations  and  nagnifi cations  nay  occur. 
Therefore,  tke  streaslines  and  equipotential  lines  of  the 
€-planc  reaain  orthogonal  in  the  plane  of  transfof ation. 

An  exception  occurs  at  points  where  =  0.  These  points 
correspond  to  flow  stagnation  points  in  the  exact  circle 
plane,  and  also  in  the  airfoil  plane,  unless  the  trailing 
edge  is  cusped,  in  which  case  the  airfoil  has  no  trailing 
edge  stagnation  point. 

Rienann  has  shown  that  the  transforned  curve  aay  be 
specified  without  negating  the  conforual  characteristic  of 
the  transformation.  In  particular,  the  interior  of  any 
sinply-connected  region  can  be  aapped  inside  a  circle,  with 
the  curve  enclosing  the  region  napping  onto  the  circunference 
of  the  circle.  This  transf ornation  is  unique  when  the  origins 
and  orientations  of  the  coordinate  systems  in  the  two  planes 
are  specified. 
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Jonkouski  Transformation 

The  transformation  £  *  z*  *  a^/z*,  where  the  C-plaie  is 
the  plane  of  .the  airfoil  (£  *  x  ♦  iy)  is  the  Joukowski  trans¬ 
formation.  If  the  airfoil  is  elliptical  or  one  of  the 
special  class  of  airfoils  known  as  Joukowski  airfoils,  then 
a  circle  results  in  the  z* -plane.  If  not,  then  a  curve 
(pseudo-circle  in  Fig.  1)  closely  approximating  a  circle 
results  for  conventional  cirfoil  shapes.  Khile  this  trans¬ 
formation  is  not  essential  to  Theodorsen's  method,  it  improves 
the  rate  of  convergence  of  the  subsequent  transformation. 

The  constant  "a"  is  included  to  preserve  dimensions  and  has 
the  dinension  of  length.  The  curve  in  the  z* -plane  can  also 
be  represented  as  z'  =  ae^+lu  in  polar  coordinates,  where 
p  =  ae^  is  the  radius  vector  and  the  angular  coordinate  is 
u.  Substituting  the  polar  form  of  z*  into  the  transformation 


equation  yields 


£  =  2a  cosb  (p  ♦  i») 


£  =  2a  cosh  p  cos  u  ♦  2ia  sinh  p  sin  ta  (3) 


Noting  that  £  =  x  ♦  iy,  and  equating  real  and  ioaginary  parts, 
the  airfoil  coorJinates  cay  be  found  explicitly  in  terns  of 
the  transforation  coordinates 


=  2a  cosh  p  cos  v 


=  2a  sinh  i)  sin  u 
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The  inverse  relations  between  coordinates  will  now  be 
determined.  Solving  Eqs  (4)  for  cosh  f  and  sinb  p  and  sub¬ 
stituting  into  the  identity  cosh2  +  -  sinh2  ♦  *  1  gives 


f  x  )  2  .  /  y  ^ 2  .  i 

\2r  cos  ct  j  ^2a  sin  «  / 


CS) 


which  can  be  solved  for  m 


=  sin"* 


(6) 


where 


a  =  1 


(^) 


Likewise,  Eqs  (4)  nay  be  solved  for  cos  «s  and  sin  w  and 
substituted  into  the  identity  cos2  w  ♦  sin2  w  =  1,  giving 


( _ ; _ ) 2 .  I 


\2a  cosh  ^2a  sinh  vy 


1*-. 


(7) 


so  that. 


#  *  s: 


inh“*  J\j2\-u  *  /n2  ♦ 


Cy/a)z) 


(8) 


Thcodorscn  Transformation  •  C„ 

__ — -  j 

n=0  z 

The  Theodorscn  transformation  z'  =  z  e  naps  the 


pseudo-circle  in  the  z-planc  into  an  exact  circle  (Fig.  1) 
in  the  z-plane.  The  coefficients  Cn  arc  complex  nuebers  of 


the  f orr.  At  ♦  ili  .  The  condition  that  the  flov  at  infiiity 


Joufcovski 
Trans formation 


jn 

tn  Thcodorsen 

Trans fornat ion 


I 


Exact-circle 


GAX/AE/72-2 


■ast  be  identical  is  airfoil  aid  exact  circle  plaaes 
requires  that  A0  *  Bc  *  0.  Tke  transformation  then  becomes 


•  C  •  A-  ♦  it_ 

I  -r  I  ~ 

a=l  *  n=l  * 


z*  *  *  e 


*  z  e 


Suitable  ueaas  cost  bow  be  fouad  for  determining  the  coa> 
staats  A)g  aad  lB. 

The  equation  for  the  exact  circle  may  be  expressed  as 
z  =  ae^*z®  ia  polar  coordiaates  where  $  is  a  coast  act  and 
K  =  ae^  is  the  radius  vector  and  the  angular  coordinate  is  6. 
Eliminating  the  constant  "a"  between  the  polar  definitions  of 
z  and  z*  yields 


(V-t)  * 


(10) 


Kov,  combining  with  the  tracs f creation  equation 


«*  A  ♦  i  B- 

(*-$)  ♦  i  (t»-e)  =  l  - - — 

n=l  z 


Note  that 


zn  =  (ae**l9)n  =  { R(cos  0  ♦  i  sin  0)Jn  = 


>s  n6  -~i  sin  n$ 


Substituting  into  F.q  fll) 
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•  A,  ♦  i*a 

(H)  *  iM)  *  I  - i -  (cos  -  i  sio  a0)  (13) 

*=1  R 


expudisj 


(♦-♦)  ♦  i(«-8) 


*  f  *»  1 

=  T  —  cos  —  sio  a0  I 

**  J 

t  r».  ^o  1 

a=l  II*  J 


(14) 


aod  equating  real  afcd  iuaginary  parts  of  Eq  (14)  establishes 
the  folloviag  conjugate  Fourier  expansions: 


-  ♦ 


-  J"  Ab  IB 

=  T  «  —  cos  n0  ♦  sin  n0 

n=l  [hn  R* 


(15) 


*  f  Bn  A* 

m  -  0  =  J  j  —  cos  n0  -  —  sin  '  9 
n=l  [R®  R 


(16) 


A  B 

The  Fourier  coefficients  —  ,  — —  ,  and  4  (the  constant  tern) 

Rn  R® 

of  Eq  (15)  are  given  by 


A  -  2r 
n  1  , 


n  "  *  cos 


A  dX 


—  /  $(A)  sin  nX  dX 

v  0 


(18) 


I?  /<,  *(*>  « 
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where  A  is  aw  angular  aeasure  corresponding  to  9  aid 
f(l)  «  ♦(••(X)J  ,  since  f  is  haova  oaly  as  a  faactioa  of  m. 
Therefore,  a  .relation  betweea  9  aad  0  (or  X)  aast  be  found. 
For  any  particular  value  of  9,  the  coefficieats  given  ia 
Eqs  (17)  aad  (li)  nay  be  substituted  iato  Eq  (16),  thereby 
eliaiaatiag  these  coefficieats  froa  that  relatioa 


I 

a=l 


2v  __ 

cos  a9  /  0(1)  sia  al  dX 


2t  _ 

sia  a0  /  0(1)  cos  aX  dX 

0 


(20) 


which  reduces  to 


u  -  6 


£  l 

n=l 


2v 


1 


/  $(X)  sia  n(X-8)dXj 

0  J 


,  _ 

i/„  ♦(*) 


J  sin  n(X-0)  dX 
n-1 


(21) 


But, 


^  sin  n(X-0) 


n=  1 


1  X-0 

J  cot  — 


cos 


[(2K+1)  - 


0  1 


2  J 


2  sin 


ITT 


(22) 


so  that 
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i 


w 


6  *  Lis 

Jf-M» 


57  /  ♦(*)  cot 

0 


1 

2* 


2*  cos  [ (2K*1) 

/  fW  - 1—  2  J 


2  sin 


X-e 


dX 


(24) 


Since  the  second  ten  is  identically  zero  and  the  first  one 
is  unaffected  in  the  linit, 

“  '  8  “  5?  /o*  *(X)  cot  dX  (25) 

The  difference,  (w-0)  or  (m-X) ,  between  angular  coordinates 
of  a  point  and  its  inage  in  the  pseudo-circle  and  exact 
circle  planes  is  called  the  confornal  angular  distortion 
function,  and  is  denoted  by  e,  so  that 

e(0)  =  /  *[e(A)  ♦  X]  cot  dX  (26) 

Equation  (26)  is  of  fundamental  importance  in  the 
practical  application  of  Theodorsen's  method,  as  it  may  be 
solved  iteratively  for  c(0),  from  which  u>(8)  is  obtained 
permitting  the  solution  of  Eq  (19)  for  the  constant  <>.  Note 
that  the  relationship  between  polar  coordinates  in  the 
pseudo-circle  and  exact  circle  planes  has  been  obtained 
without  having  to  solve  for  the  Fourier  coefficients  directly. 
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Surface  Velocity 

It  is  shown  in  potential  flow  theory  that  the  complex 
velocity  potential  for  a  circle  in  the  z-plane  immersed  in 
steady,  x-directed  flow  is 

w(z)  *  $  ♦  iV  =  -V^  ^  ♦  f“)  -  T?  In  I  (27) 

where  V„,  is  the  freestream  velocity,  R  =  ae^  is  the  radius 
of  the  circle,  and  T  is  the  circulation.  The  fluid  velocity 
in  the  circle  plane  is  obtained  by  differentiating  w(z) 

V(z)  •  HI  =  u  +  iv  B  “Vw  ^1  -  pr)  -  H  (28) 

This  velocity  will  be  zero  (u  =  v  *  0)  at  stagnation  points. 

Therefore,  the  locations  (zc)  of  these  points  is  given  by  the 

value  for  z  when  ^  s  fl.  That  is, 

dz 


A  non-zero  angle  of  attack  a  may  be  accounted  for  by 
replacing  z  by  zeict  in  the  real  part  of  Eq  (27)  .  This  has 
the  effect  of  rotating  the  flow  field  about  the  circle  by 
the  angle  a. 

w(z)  =  -V„  |  zela+  e"ia>J  -  In  z  (30) 
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O 


and 

r  V.(I)  *  -v„  eia  (  1  -  si  e-2i«  j  -  il_  (3D 

The  well-known  Kutta- Joukowski  condition  requires  that 
the  trailing  edge  (w  =  it)  on  an  airfoil  which  does  not  have 
a  cusped  trailing  edge  be  a  stagnation  point  in  order  that 
infinite  velocities  and  pressure  gradients  may  be  avoided. 

The  circulation  required  to  meet  this  condition  may  be  found 
by  setting  z0  *  Re*^Tr+^  °  -Re*^  (where  8  is  the  negative  of 
the  value  of  e  at  u  ■  f,  and  is  called  the  angle  of  zero  lift 
for  the  airfoil)  and  V(z)  *  0  in  Eq  (31),  and  solving  for 
the  corresponding  T: 

T  =  4ttRVco  sin  (a+B)  (32) 

Then 

V(z)  .  pL  .  -V„ei0  (l  -  4  *-aa  )  -  l  41,RV"  Sln.  (a<6)  (33) 
dz  \  z*  I  2ttz 


The  velocity,  v(z),  in  the  plane  of  the  circle  must  now 
be  converted  to  the  velocity,  v(C)  ,  at  the  surface  of  the 
airfoil  through  the  use  of  transformation  derivatives: 


v  (C) 


dw  d_z  dz 1 
dz  dz'  d£ 


(34) 


and  these  derivatives  must  now  be  determined.  Recall  that 
z  =  Re*®  on  the  surface  of  the  circle.  Substituting  this 
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value  for  z  in  Eq  (33)  gives 

37  *  elaJl  -  e'2l(c*e))  -  2iV— e-i®  sin  (a*£) 
which  reduces  to 

=  -2iV„e'ie[sin  («♦«)  ♦  sin  (o*6)J 
d  z 

the  absolute  value  of  which  is 

j^|  =  2Y„[sin  (b+8)  *  sin  (a*&)) 

Fron  the  definitions  of  z,  z ' ,  and  C 

z*  _  ae^**”  _  (tf-0)*ie 

z  =  e+i6 
ae  ‘ 

Differentiating  and  noting  that  $  is  a  constant 

Ifl  =  ^  ((S-«*icl  = 

r-  u  4  z  ji  c**ioi 

But 

1  dz  =  dlae^*10)  =  idO 

2  ac**10 

Substitution  into  I:q  (39)  gives 


(39) 


(40) 
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dz* 

JT" 


ac^ic) 

T3d 


(41) 


Further  algebraic  aanipalatioa  oa  this  equatioa  gi res 


dz*  z* 

3T"  =  — 


1 

cae/a««)  ♦  (ac/a**)  -  i  (a*/a«)  i 

- (W/dS) - 


(42) 


Appljing  the  definition  c  *  w  -  0,  aad  reducing 


az* 

az 


ll 

z 


1  -  i  (d$/du) 

1  -  ( dc/dw) 


1 

[ 

I 

J 


(4S) 


And  finding  the  absolute  value 


IJilL  .*-* 

/l  ♦  (dy/dts)2  ] 

Idz  |-  e 

1  -  (dc/dw)  j 

(44) 


Inverting 


dz 

dz' 


=  e*'* 


fde/dfc)) 


/TT 


ca^/d«)z 


(45) 


which  is  the  desired  derivative.  Substituting  the  definition 
z'  *  ae^+*°  into  the  Joukowski  transfornation 


C  =  ae^+iu  ♦  ae”^_lu  =  2a  cosh  (*+iw)  (46) 


Differentiating 


dr._ 

dz  • 


2  si nh 


(V»+i«i) 


e-*“iw 


(47) 
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Expaa4iif 

3§r  -  2e"^"i<#  (siah  ♦  cosh  i-  ♦  cosh  +  siah  i«)  (4t) 

■otiaf  that  cosh  in  *  cos  «  aad  siah  ii*  i  sia  m,  this 
hecoaes 


||r  *  2e”*'i-  (siah  ♦  cos  a  *  i  cosh  f  sia  «)  (49) 

Thea 

Jgij.|  *  4e'2^  (siah2  f  cos 2  a  ♦  cosh2  f  sia2  m)  (50) 
which  reduces  to 

|g|r|  *  4e~2*  (si ch2  ♦  ♦  sia2  a)  (51) 

Taking  the  square  root  and  inverting  gives  the  final 
derivative 


dz 1 
dC 


e* 


2/sinh^  $ 


♦  sm*  v 


(52) 


The  derivatives  given  in  Eqs  (37),  (45),  and  (52)  will 
now  be  substituted  into  Eq  (34)  to  find  the  airfoil  surface 
velocity 


V  e^fsin  (a+8)  ♦  sin  (a*B)][l  -  (de/du)] 

V(0  B  -  - - - 

/(sinh2  ♦  sin2  w)[l  ♦  (d^/dto) 2 ) 


(S3) 
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III.  Kanttan-Fohlhausen  Method 


Origin 

Geaeral  aotioa  ef  a  viscous  fluid  is  described  hjr  the 
■avier-Stokes  equations,  which  were  derived  froa  first 
principles.  Schlichting  (lef  3)  has  showa  that,  if  Keynolds 
aaaber  is  large,  viscosity  is  coastaat,  aad  flow  is  steady 
aad  iacoapressible,  thea  the  eqnatioas  of  aotioa  for  flow 
over  a  surface  aay  be  approximated  by  the  aoaeatua  iategral 
equal ioa: 


t  All  Tw 

u  as-  *  <2*2  +  *0  u  air  -  pr  <54> 

where  the  boundary  layer  displaceaent  and  aonentua  thick¬ 
nesses  are,  respectively. 


(55) 


(56) 


and  where  x  and  y  are,  respectively,  the  streanwise  and 
cross-stream  coordinates  of  a  general  point  in  the  flow 
(Fig.  2),  with  the  corresponding  velocities  given  by  u  and 
v.  Note  that  the  boundary  layer  edge  velocity  U  is  a 
parameter  common  to  both  the  potential  flow  and  boundary 
layer  regions. 

The  following  analysis  is  due  originally  to  Pohlhausen, 
but  is  presented  in  the  form  developed  by  Holstein  and 
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Mira,  aoA  described  by  Scklicktiif. 

TeMity  Profile 

A  velocity  profile  of  the  following  font  is  ion  asme4: 

}  •  an  «  M2  ♦  en3  ♦  4n4  (57) 

where  n  *  jr/4  aai  a,  h,  c,  aod  A  are  coefficieats  to  he 
AeteroiaeA . 

Observatioas  iaAicate  that  the  followiag  boundary  coodi- 
tioas  shoalA  he  applied  to  the  solntioa  of  Eq  (54): 


a  *  0 

at  y  -  0 

(SS) 

at  y  *  0 

(59) 

o*0 

at  y  *  6 

(60) 

at  y  =  6 

(61) 

i-S  -  o 

>y 

at  y  -  6 

(62) 

where  6,  called  the  boundary  layer  thickness,  is  the  cross- 
strean  distance  froa  the  surface  at  which  the  boundary  layer 
aeets  the  potential  flow  region. 

Note  that  the  first  boundary  condition,  Eq  (S8) ,  is 
identically  satisfied  by  the  assuecd  solution  by  virtue  of 
the  omission  cf  a  constant  tern  in  the  polynomial.  There 
remain,  therefore,  four  boundary  conditions  to  determine 
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a,  h.  c.  a»i  4.  I»  tens  af  the  assaacrf  velocity  profile, 
tic  sacaaJ  thravfl  fifth  hevadary  ceeditieas  heceoe. 
respective  ly,- 


at  n  ,  •  (63) 

at  a  *  1  (M) 

at  «i  =  1  (65) 

at  v>  =  1  (ii) 


Applicatioa  of  these  coed  it  ions  to  the  asstsned  solntioa  aad 
defining  a  shape  factor 


6^  dU 
v  dx 


(67) 


leads  to  the  following  expression  for  the  velocity  profile: 
£  =  (2n  -  2n3  ♦  n4)  ♦  ^  (n  -  3n2  ♦  3n^  -  n4)  (68) 

U  6 

3u 

By  applying  the  physical  constraints  that  jy  =  0  at 
y  =  o  at  the  separation  point  and  that  ^  <  1  within  the 
boundary  layer,  the  upper  and  lower  physical  Units  on  A 
arc  readily  deternined.  These  have  been  found  to  be  *12 
and  -  12 ,  respectively. 


CAJf/AE/72-2 


Other  Parameters 

The  displacement  thickness,  momenta*  thickness,  and 
shear  stress  mill  non  he  foand  in  terns  of  A. 

Substituting  the  a  bore  velocity  prefile  and  the 
definition  for  n  into  the  definition  for  displacement  thick¬ 
ness  given  by  tq  (SS)  gives 

A.  , 

X~  *  /  [1  -  (2ri  -  2rr3  ♦  *,<)  - 

•  0 

|  («i  -  3n2  ♦  3n3  -  n^Jldn  (*9) 


And.  finally. 


61  m  3_ 

«T  10  120 


(70) 


Likewise. 

*2  _  37  A  _  A2 

3T"  3l?  945  9072 

Kov.  the  shear  stress,  by  definition,  is 


(71) 


T 


w  _ 


y 


(72) 


Evaluating 

gives 


the  derivative  and  nultiplying  both  sides  by 


6 

yu 


T  ft 


A 
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Siaplificatioa  *f 


Iitetnl  ImatiM 


The  MwatM  iategral  eqwatioa  caa  aw  he  pat  late  a 
fora  fran  which  practical  resales  aajr  he  obtained. 

Naltip lying  hath  sides  ef  Eq  (54)  by  ,  aad  redaciag. 


leads  te 


tw*2  V  .  **2  /  *1  \*2*  dB 

~W  *  a  *2  J5~  ♦  dx 


(74) 


Eqaatioas  (67),  (6t),  (76),  (71),  aad  (74),  together 
with  the  haewa  edge  velocity  distribatioa,  fon  a  solvable 


syst« 


However,  eliwieatioa  of  deoeadeace  of  the  other 


variables  oa  i  is  desirable,  siace  the  bowadary  layer  thi.ck- 

d6y 

aess  is  physically  aeaaiagl ess-  Fortheraore,  the 

quaatity  ia  Eq  (74)  leads  to  a  reqeireaeat  for 
d^U 

— j  at  every  peiat  oa  the  sarface.  By  sai table  aatbeaatical 

dX  d2  U 

aaaipulatioa,  this  caa  be  reduced  to  a  reqoireaeat  for  -~ 

dr2 

at  the  upstreaw  stagflation  poiat  only.  The  following  analysis 
proceeds,  therefore,  toward  these  objectives. 

Introdoction  of  several  aew  paraaet ers  will  perait 
siaplificatioa  of  the  fora  of  Eq  (74).  Define 


i  2 

*2  dU 
V  dx 


(7S) 


Then , 


r  -  I6?  \  2  62  i  (  62  \  2  .  ("SI 

t-\r) 


_A_  -  A__ 

945  9072 


(76) 
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So  that 


2  ,  d62 

vJ2dT 


(78) 


Also,  from  Eqs  (75)  arid  (77) 


K  =  Z 


(79) 


The  ratio  of  displaceaent  to  non en tun  boundary  layer  thick- 
nesses  nay  be  defined 


f  -  1 

1 


(«0) 


Or,  using  Eqs  (70)  and  (71) 


3  A 
10  "  120 

37  t  A2 
315  "  945  '  9072 


(81) 


And,  finally,  the  definition 


f  =  lifi 

*2  HIT 


(82) 


gives,  using  Oqs  (71)  and  (73), 


77 


X  /« 

f,  *  — £ 

2  ||U 


(/)  -  (2  *  2)  ( sis  - 


A  A2  . .... 
945  “  5072  I  l83J 


) 


Mow,  if  these  definitions  are  substituted  into  Eq  (74) , 
the  result  is 


U  dZ 


2  «  J  £  +  <2  *  fl)K 


(84) 


which  aay  be  rearranged  to  give 


"  E  *  2f2  -  «  ♦  «l)* 


(*S) 


The  right  side  of  Eq  (8S)  is  a  function  of  A  alone,  as  shown 
by  Eqs  (76),  (81),  and  (83).  Therefore,  if  a  new  parameter 
F  is  defined  such  that 


F  =  2f2  -  (4  ♦  2fj)K 


(86) 


Then,  it  nay  be  shown  that 


F  c  2  I57  -  A  -  *2  \ 

h  ^3T5T  945  9072  / 


116 

315 


A  ♦  (  2  ♦  2  JA2  ♦ 

\ 945  120  ) 


2  A3 


9072 


(87) 


And  Eq  (85),  the  equation  of  motion,  becomes 


dZ  F 
HI  =  U 


(88) 
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Conditions  at  Stagnation  Point 

Enough  relations  now  exist  to  permit  a  stepwise  solution 
for  all  the  boundary  layer  pararaeteis,  once  a  starting  condi¬ 
tion  is  known.  Starting  conditions,  therefore,  will  now  be 
determined  at  the  upstream  stagnation  point. 

Since  the  flow  is  brought  to  rest  at  the  stagnation 
point,  Uc  s  0  (where  the  "o'*  subscript  means  that  an  asso¬ 
ciated  variable  is  to  be  evaluated  at  stagnation  point 

conditions).  Then,  from  Eq  (88),  F0  must  also  be  zero  if 
(  dZ  \  I 

^  <Tx  j  |0  is  to  exist.  Therefore,  AQ  must  be  a  root  of 
Eq  (87)  when  F  =*  0 .  The  roots  of  this  equation  are  -72.2, 
-37.8,  7.052,  17.9,  and  28.2.  Because  of  the  physical 
limitations  on  the  range  of  A,  it  may  be  seen  that 


AQ  =  7.052 


Then  from  Eq  (76) , 


.0770 


And,  from  Eq  (79), 


Zo  a 


.0770 

(S)l. 


which  determines  ZQ,  since  the  edge  velocity  distribution  is 
known  from  the  potential  flow  analysis. 

Now,  (  3Y  )  may  f°unt*  by  applying  L'Hospital's  Rule 
to  Eq  (88) 
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(8)1. ' 


t4_  (dF/dx)  (dF/dKHdK/dx) 

(dU/dx)  =  (dU7dxJ 


(92) 


But,  from  Eq  (79) , 


dK  _  dZ  dl)  .  7  d2U 

HIT  33T  53T  4  ^2 


Substituting  this  expression  into  Eq  (92)  and  passing  to  the 


limit  gives 


)  lo  ■ 
ai'ii  (<un 

)  lo  l  a*  /  o 


But,  using  Eqs  (76)  and  (87), 


(ay)  L 

(ay)lo 


=  -5.57 


Substituting  the  values  for  ZQ  and  given  by  Eqs  (91) 

and  (95),  respectively,  into  Eq  (94)  provides  an  expression 
for  in  terms  of  the  edge  velocity  distribution: 


=  -.0652 
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Bsiijt  the  Method 

With  starting  conditions  known  at  the  stagnation  point, 
the  value  for-  Z  at  a  distance  Ax  dcvnstrean  is 


Z,  = 


7  /  <*2  \ 

z°*  l3*/ 


Ax 


(97) 


Then,  from  Eq  (79), 


(98) 


The  corresponding  value  of  A,  fron  Eq  (76)  may  now  be  used 
to  find  Fj  fron  Eq  (87).  Using  Eq  (88), 


dZ  ) 

**/  1 


1 


u, 


(99) 


so  that  the  stepping  process  nay  be  continued  indefinitely, 
the  value  of  Z  at  point  "n"  being 


Zn  =  Zn- 1  * 


dZ 

dx 


Ax 


(100) 


n- 1 


Knowledge  of  An  permits  determination  of  fj^,  ^2n»  and  t^e 


u  \ 


,  using  Eqs  (45),  (83),  and  (68), 


velocity  profile  j 
respectively.  Substituting  Zn  into  Eq  (77)  yields  ^2n» 
which  in  turn  produces  from  Eq  (82),  and,  the  ultimate 

objective,  6^  from  Eq  (80).  This  value  of  must  be 

added  to  the  original  airfoil  thickness  to  form  a  better 
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approxinati on  of  the  potential  flow  airfoil  shape  to  be 
analyzed  again  by  Theodor sen's  Method. 

Linitation 

The  Karnan-Pohlhausen  nethod  of  boundary  layer  analysis 
does  not  properly  describe  the  flow  downstrean  of  a  separa¬ 
tion  point  (the  point  at  which  the  shear  stress  at  the  wall 
vanishes) . 
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IV.  Conputcr  Study 
Purpose  of  Coaputer  Stodv 

A  Fortran  Extended  progran  for  the  CDC-6600  con pater 
was  written  as  part  of  this  study  and  is  listed  in  Appendix 
A.  It  inplenents  the  theory  discussed  in  the  preceding 
chapters  and  lays  a  foundation  for  nore  sophisticated  air¬ 
foil  analysis  prograns.  Flow  paraaeters  fron  this  progran 
nay  be  conpared  with  those  fron  sinplified  airfoil  theories 
in  order  to  deteinine  the  validity  of  thin  airfoil  and  in- 
viscid  assumptions.  Sone  pertinent  facts  about  the  progran 
will  be  discussed  in  this  chapter. 

Progran  Composition 

The  computer  program  is  composed  of  a  main  progran  and 
several  subprograms  and  functions.  The  main  progran  is 
called  MAGIC  and  provides  the  means  for  reading  input  data, 
controlling  data  flow  and  sequencing  among  subprograms,  and 
testing  for  completion  of  the  iteration  on  displacement 
thickness.  One  principal  subprogram  is  THEO,  which  applies 
Theodorsen's  method  to  the  airfoil;  the  other  is  BOUND, 
which  solves  the  boundary  layer  equations  using  the  Karman- 
Pohlhausen  method.  Other  subprograms  and  their  purposes 
include:  FNEVAL,  which  defines  integrands;  SIMPS,  which 

performs  integrations  using  Simpson's  Rule;  MTXEQ,  which 
solves  systems  of  alpcbraic  equations;  PLSQ,  a  polynomial 
least-square  curve-fitting  routine;  ATKN ,  an  interpolating 
function;  several  subroutines  for  the  CALCOMP  plotter;  and 
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uaerovs  standard  functions  frou  the  system  library. 

Input 

Pro gran  MAGIC  requires  an  alphanumeric  description  of 
the  airfoil  (e .  g . ,  NACA  1408),  the  number  of  data  points, 
the  Cartesian  coordinates  of  each  point,  a  plotting  index, 
freestream  velocity,  viscosity,  density,  and  the  angle  of 
attack  on  the  chord  line.  The  airfoil  chord  line  is  assumed 
to  lie  on  the  x-axis  with  the  trailing  edge  in  the  positive 
x-direction  and  the  midchord  point  at  the  origin.  Joints 
defining  the  airfoil  surface  must  begin  at  the  trailing  edge 
and  be  nuubered  in  the  counterclockwise  direction.  All 
dimensions  must  be  referenced  to  chord  length.  The  plotting 
index  causes  the  airfoil  and  its  transformations  in  sub¬ 
program  THEO:  (1)  to  be  plotted  on  the  same  axes,  (2)  not 
to  be  plotted,  or  (3)  to  be  plotted  on  separate  axes, 
according  as  its  value  is  -1,  0,  or  +1,  respectively.  The 
angle  of  attack  nust  be  given  in  degrees. 

Additionally,  accuracy  requirements  on  various  computa¬ 
tions  may  be  reset  by  modifications  to  the  program. 

Output 

Program  MAGIC  computes  and  prints  Cartesian  coordinates 
in  all  three  planes,  polar  coordinates  in  the  pseudo-circle 
and  exact  circle  planes,  potential  flow  surface  velocity, 
coefficients  for  fifth-order  polynomials  for  ^(6)  and  $(u))  , 
radius  of  exact  circle,  zero-lift,  angle  of  attack,  circula¬ 
tion,  section  lift,  pressure  coefficients,  shear  stress  at 
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surface,  ■neitra  thickness,  displacement  thickness,  shape 
factor,  aad  several  iaportaat  quantities  ased  ia  computing 
these.  The  airfoil  surface  ia  all  three  plaaes  is  plotted 
(provided  that  the  plotting  iadex  is  aoa-zero) ,  as  are  pres- 
sore  coefficient,  shear  stress,  and  displaceaeat  thickaess 
as  fnactions  of  chordwise  statioa.  Saaples  of  ootpot  format 
are  shown  ia  Appendix  A. 

Limitations  on  Program 

Although  subprogram  THEO  gives  correct  results  for  non¬ 
zero  angles  of  attack,  subprogram  BOUND  does  not,  since  the 
shift  of  upstream  stagnation  point  with  angle  of  attack  has 
not  been  accounted  for.  Therefore,  accurate  results  for  the 
complete  program  nay  be  expected  only  for  zero  angle  of 
attack. 

Values  of  shear  stress  and  displacement  thickness 
generated  by  subprogram  BOUND  downstream  of  the  separation 
point  becone  very  large.  These  values  are  printed  as  output, 
but  are  actually  United  to  a  reasonable  naximun  value 
before  being  plotted  or  transmitted  to  subprogram  THEO  for 
another  iteration.  In  particular,  the  displacement  thick¬ 
ness  at  the  trailing  edge  is  set  to  zero  in  order  that  the 
new  airfoil  shape  considered  by  subprogram  THEO  will  be  a 
closed  curve. 
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X.  Results 

Couforual  Transformations 

Tramsforuat ions  were  accoaplisheJ  by  sskprofraa  TBEO 
oa  several  airfoils,  including  XACA  1408  and  XACA  4424. 

Tie  output  plot  of  the  traasfomatioa  oa  the  KACA  4424 
airfoil  is  presented  ia  Fig.  3.  Success  of  the  transforma- 
tioa  oa  these  airfoils  was  ueasured  by  the  accuracy  of  the 
output  angle  of  zero  lift,  £ .  The  values  for  these  angles 
as  given  by  Abbott  2nd  von  Doenhoff  (Kef  1)  are  presented 
for  conparison  in  Table  1.  Conparison  data  was  not  available 
for  the  other  airfoils  used. 

Table  1 

Conpzri son  of  Zero  Lift  Angles  of  Attack 

Airfoi 1  (Progran  MAGIC)  (Ref  S) 

KACA  1408  .95*  .9S# 

KACA  4424  3.58*  3.50* 

Potential  Flow  Surface  Velocities 

Surface  velocities  generated  by  subprogran  THEO  were 
also  compared  with  those  given  by  Abbott  and  von  Doenhoff. 

A  typical  conparison  is  shov;n  in  Table  2.  The  velocities 
agree  within  2"  in  most  cases.  Exceptions  occur  at  points 
where  the  radius  of  curvature  is  relatively  snail;  that  is, 
near  the  leading  edge.  This  probably  reflects  poor 
accuracies  in  nurcrical  derivatives  in  that  area. 
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Table  2 

Ccftriswi  >f  Poteatial  flaw  Velocities 
for  X.ACA  1401  Airfoil 

Upper  Sarface  Lower  Sor face 


Chordwise 

Statioa 

(perceat 

chord) 

(Prograa 

MAGIC) 

Y/Ym 

C*ef  1) 

Chord vise 
Statioa 
(Perceat 
chord) 

T/Tm 

(Prograa 

MAGIC) 

v/*m 
(*ef  1) 

0 

0 

0 

0 

0 

0 

1.189 

.780 

1.061 

1.311 

.829 

1.039 

2.413 

1.067 

1.121 

2.582 

1.079 

1.0S9 

4.896 

1.135 

1.151 

S.  104 

1.121 

1.105 

7.3S6 

1.142 

1.161 

7.614 

1.112 

1.105 

9.883 

1.145 

1.161 

10.117 

1.102 

1.099 

14.885 

1.153 

1.164 

15.111 

1.096 

1.092 

19.904 

1 . 152 

1.162 

20.096 

1..085 

1.082 

24.926 

1.148 

1.157 

25.074 

1.075 

1.071 

29.950 

1.139 

1.149 

30.050 

1.065 

1.063 

40.000 

1.126 

1.131 

40.000 

1.054 

1.047 

SO. 020 

1.107 

1.110 

49.980 

1.038 

1.034 

60.034 

1.087 

1.089 

59.966 

1.023 

1.019 

70.041 

1.067 

1.070 

69.959 

1.012 

1.008 

80.039 

1.042 

1.043 

79.961 

.998 

.991 

90.027 

1.014 

1.003 

89.973 

.982 

.965 

95.016 

.  966 

.983 

94.984 

.953 

.955 

100.000 

0 

0 

100.000 

0 

0 
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lowliry  Layer  Parameters 

Pressure  coefficients,  shear  stress,  and  boundary  layer 
displacement  thickness  distributions  are  qualitatively  correct 
(Figs.  4,  5,  and  6)  upstreaa  of  the  separation  point.  As  a 
quantitative  check,  the  separation  point  for  uniform  flow 
around  a  circle  was  computed.  Subprogram  BOUND  found  this 
point  to  be  at  105*  from  the  upstream  stagnation  point.  This 
agrees  well  with  values  from  104.5*  to  108.8*  determined  by 
various  investigators,  as  cited  by  Schlichting. 

Running  Time  and  Core  Memory  Required 

Computer  tine  required  for  program  MAGIC  is  a  function 
of  airfoil  shape,  number  of  points,  and  accuracy  required. 

The  KACA  1408  airfoil,  with  34  points  defined,  required 
approximately  50  seconds  of  CDC-6600  central  processor  time 
and  two  iterations  to  obtain  a  displacement  thickness 
accuracy  of  .1%  of  chord  length. 

Approximately  40,000  octal  words  of  core  memory  were 


required. 
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VI.  Conclusions  and  Recommendations 

Conclusions 

Based  upon  the  results  shown  in  the  preceding  chapter, 
the  following  conclusions  were  drawn  from  this  study: 

1.  Theodorsen's  method  computes  values  of  zero  lift 
angle  of  attack  and  potential  flow  surface  velocities  very 
accurately. 

2.  The  Karman-Pohlhausen  method  of  solving  the 
boundary  layer  equations  yields  good  results  upstream  of  the 
separation  point.  No  useful  information  is  obtained  from 
these  equations  downstream  of  the  separation  point. 

Knowledge  of  the  shape  factor  permits  determination  of  the 
coefficients  in  the  assumed  velocity  profile  at  any  point 
along  the  surface. 

3.  Program  MAGIC  is  a  practical  means  of  analyzing 
incompressible  flow  about  an  airfoil,  without  resort  to  small 
perturbation  or  inviscid  assumptions.  Core  memory  and 
central  processor  requirements  are  very  reasonable  for  an 
analysis  of  this  magnitude.  Input  data  requirements  are 
minimal . 

Recommendations 

The  following  recommendations  are  made  for  future 
studies  in  extension  of  this  one: 

1.  Program  BOUND  could  be  modified  to  permit  boundary 
layer  computations  when  angle  of  attack  is  not  zero. 
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2.  Integration  of  pressure  distribution  and  shear 
stress  would  provide  drag,  moment,  and  center  of  pressure 
information. 

3.  According  to  Giesing  (Ref  2),  Theodorsen's  method 
can  be  extended  to  linear  arrangements  of  identical  airfoils 
Program  MAGIC,  so  extended,  could  be  of  considerable  value 
in  cascade  theory. 

4.  Theodorsen  states  that  accounting  for  leading  edge 
and  trailing  edge  radii  in  locating  the  coordinate  system 
origin  improves  the  rate  of  convergence  of  the  iteration 

on  e .  Inclusion  of  this  correction  might  result  in  some 
saving  in  computer  time. 

5.  Numerical  derivatives  of  a  higher  order  of  accuracy 
near  the  leading  edge  would  probably  improve  velocities 
computed  in  that  region. 

6.  The  inverse  problem  (see  Ref  4)  could  be  programmed 
That  is,  given  a  desired  pressure  distribution,  the  required 
airfoil  shape  could  be  determined. 

7.  A  more  graphic  means  of  input  and  output  could  be 
obtained  by  adapting  the  program  to  the  cathode  ray  tube 
display  terminal. 
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Appendix  A 


Computer  Program  and  Sample 


Format 


Reproduced  from 
best  available  copy. 
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Ellie  B.  Underwood,  Jr.,  was  born  on  27  January  1939 
in  tieevil le,  Texas.  He  was  graduated  in  June  1956  froa 
Ballinger  High  School,  Ballinger,  Texas.  In  June  1961,  he 
was  awarded  a  Bachelor  of  Science  degree  in  Aerospace 
Engineering  froa  the  University  of  Texas,  Austin,  Texas,  and 
was  coaaissioned  in  the  United  States  Air  Force  at  the  saae 
tine.  In  June  1960,  he  was  assigned  to  the  Air  Force 
Institute  of  Technology,  Wright-PatteTSon  Air  Force  Base, 
Ohio. 
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